First-principles calculations of structure and elasticity of hydrous fayalite under high pressure
Zhang Chuan-Yu1, †, Wang Xu-Ben1, Zhao Xiao-Feng1, Chen Xing-Run1, Yu You2, Tian Xiao-Feng3
College of Geophysics, Chengdu University of Technology, Chengdu 610059, China
College of Optoelectronic Technology, Chengdu University of Information Technology, Chengdu 610225, China
College of Nuclear Technology and Automation Engineering, Chengdu University of Technology, Chengdu 610059, China

 

† Corresponding author. E-mail: zhangchuanyu10@cdut.cn

Project supported by the National Natural Science Foundation of China (Grant Nos. 11404042 and 11604029), the Research Fund for the Doctoral Program of Higher Education of China (Grant No. 20135122120010), and the Open Research Fund of Computational Physics Key Laboratory of Sichuan Province, Yibin University (Grant No. JSWL2015KFZ02).

Abstract

The structures, elasticities, sound velocities, and electronic properties of anhydrous and hydrous fayalite (Fe2SiO4 and Fe1.75H0.5SiO4) under high pressure have been investigated by means of the density functional theory within the generalized gradient approximation (GGA) with the on-site Coulomb energy being taken into account (GGA+U). The optimized results show that H atoms prefer to substitute Fe atoms in the Fe1 site. Compared with the anhydrous fayalite Fe2SiO4, the mass density, elastic moduli, and sound velocities of Fe1.75H0.5SiO4 slightly decrease. According to our data, adding 2.3 wt% water into fayalite leads to reductions of compressional and shear wave velocities ( and ) by 3.4%–7.5% and 0.3%–3.4% at pressures from 0 GPa to 25 GPa, respectively, which are basically in agreement with the 2%–5% reductions of sound velocity obtained by the experimental measurement in the low velocity zones (LVZ). Based on the electronic structure, the valence and conduction bands are slightly broader for hydrous fayalite. However, hydrous fayalite keeps the insulation characteristics under the pressures up to 30 GPa, which indicates that hydration has little effect on its electronic structure.

1. Introduction

It is well known that water in the earth’s mantle exists in hydrated metamorphic minerals, anhydrous silicate minerals and dense hydrous magnesium silicate (DHMS). The importance of nominally anhydrous minerals (NAMs) incorporating a little hydrogen as lattice defects in the mantle was reported forty years ago.[1] Water of mantle minerals exists in the form of hydroxyl (OH), which had been confirmed in olivine.[2] In recent years, many researchers were concerned about the important effects of hydration of NAMS on geophysical processes and geochemical processes.[311]

The high-pressure behavior of olivine (MgxFe1−x)2SiO4 is the focus of attention, which is helpful to understand the compositions, structure, and dynamics of deep earth. In particular, olivine could be used to explore the origin of the discontinuities of density, seismic wave velocities, and electric conductivity observed at 410-km and 660-km depth inside the earth.[1217] Many experimental results on the hydrous properties of olivine have been reported in recent years. In particular, electric conductivity, elastic property, and sound velocity were significantly affected by the presence of water.[1823] It is known that water could reduce the elastic moduli of olivine under high pressure conditions.[2426] However, the elastic properties of hydrous olivine had been rarely reported. For anhydrous fayalite, the elastic moduli of fayalite have been measured by means of static pressure method.[2730] The elastic properties of anhydrous fayalite have been obtained under both high pressure and high temperature conditions.[3134] The wave velocities of fayalite in olivine-spinel phase transition were obtained employing x-ray diffraction technique.[35]

Recently, quantum mechanics method based on the first principle calculation has become an important technology to investigate the mineral material properties under extreme conditions.[3643] However, the influences of hydration on elasticities and electronic structures of olivine were rarely reported through atomistic modeling. Only structure, elastic, and electronic properties of hydrous forsterite under high pressure were explored by first-principle simulation.[4447] The density of state of electrons in hydrous forsterite was calculated using the density functional approach.[48] In addition, the electronic structures, phase transition, elastic properties, and magnetism of anhydrous fayalite had been studied employing the first-principle modeling approach.[4951] The elasticities of (Mg0.875Fe0.125)2SiO4 crystals under high pressure conditions were also calculated.[52]

Theoretical data of elastic constants of hydrous fayalite under high-pressure conditions have rarely been reported so far. The elasticity and sound velocity of hydrous fayalite could be used to assess the characteristics of seismic wave anomaly in the mantle. In this article, the hydration effect on the electronic structures, elasticities and sound velocities of fayalite under pressures up to 30 GPa has been investigated in detail.

2. Computational details
2.1. Computational method

All calculations of Fe2SiO4 and Fe1.75H0.5SiO4 were based on the first-principles method employing the VASP program.[53] To solve the electron exchange and correlation correction, the generalized gradient approximation (GGA) of Perdew–Burke–Ernzerhof scheme is introduced.[54] Due to 3d elections of Fe atom, our calculation incorporates a simplified GGA+U scheme (U = 4.5 eV and J = 0.90 eV for Fe2+). According to the test data, the cutoff energy and k-points of the irreducible wedge of the first BZ is 500 eV and 4×2×3, respectively. The total energy and the force components on each atom are converged to within 1.0×10−5 eV/atom and 1.0×10−3 eV/Å, respectively.

2.2. Computational models

Fe2SiO4 is described in the orthorhombic lattice structure with Pbnm space group. According to the previous data, H atom could occupy the cationic vacancy due to the lattice defect. To design the model of hydrous fayalite, Fe2SiO4 is transformed into Fe1.75H0.5SiO4 by replacing an Fe atom by two H atoms in the unit cell (Z = 4, 28 atoms). There are two non-equivalent sites for Fe atoms in fayalite crystal lattice, which are marked as Fe1 and Fe2, respectively. Therefore, H atoms can replace Fe atoms at two different positions. In this paper, we design six possible models for the hydrous fayalite (Fe1.75H0.5SiO4), through which an H atom and one of the O atoms forms an OH group outside the SiO4 tetrahedron. Based on the stable structure of anhydrous fayalite under the 20 GPa pressure, five proposed structures (M1, M2, M3, M4, and M5) lie in two H atoms occupy Fe1 site, while M structure catch out to be Fe2’s occupation by the H atoms, as shown in the following Fig. 1.

Fig. 1. (color online) Crystal structures of Six calculated models for hydrous fayalite Fe1.75H0.5SiO4.
3. Results and discussions
3.1. Crystal structures of Fe2SiO4 and Fe1.75H0.5SiO4

To evaluate the accuracy of the first-principles method, the crystal parameters of anhydrous fayalite Fe2SiO4 were firstly optimized. Theoretical and experimental data of Fe2SiO4 crystal are shown in Table 1. It is obvious that our calculated data are in good agreement with those theoretical and experimental results.[29,39,50,51,55] At 25.4 GPa, the value deviates by less than 1.3%.[56] The discrepancy in cell volume of the calculated and experimental values is only about 2.5%, and the maximum deviation of lattice constants is less than 0.08 Å (Δa = 0.050 Å, Δb = 0.011 Å, Δc = 0.076 Å). To compensate for the temperature induced errors, we calculated the lattice parameters and mass density of Fe2SiO4 under 3 GPa pressure condition, as shown in Table 1. It indicates that the 3 GPa pressure correction brings good agreement between computation and experiment.[29,55] Compared with the experimental data, the absolute shifts of lattice constants are less than 0.04 Å (Δa = 0.031 Å, Δb = 0.04 Å, Δc = 0.038 Å), and the difference of mass density is only 0.01 g/cm3. It is obvious that first principles calculation is an accurate method, which could be extended to calculate the properties of Fe1.75H0.5SiO4.

Table 1.

The experimental and calculated lattice parameters and mass densities of Fe2SiO4 under different pressures.

.

To obtain the stable crystal structure, six designed models of Fe1.75H0.5SiO4 at 20 GPa have been optimized. The experimental and theoretical investigations of hydrous fayalite have been rarely reported so far, thus the crystal symmetry is set to P1 in the whole simulation process. All these optimized data are listed in Table 2. The total energy of M1 model is lower than that of the other five models, which indicate that M1 model is the most stable structure among all these six models. It means that H atoms prefer to occur at the Fe1 site for the hydrous fayalite, which is consistent with the experimental results of hydrous forsterite.[57] The calculated data of anhydrous fayalite at 20 GPa are: a = 4.753 Å, b = 9.926 Å, c = 5.941 Å, and V = 280.31 Å3. For Fe1.75H0.5SiO4, the large decrease (0.081 Å) in the a axis and the relatively small decrease (0.058 Å) in the c axis are obtained, while the cell parameter b value increases by 0.5% in contrast with Fe2SiO4. The overall lattice compression of Fe1.75H0.5SiO4 is about 2.3%, which is due to the difference atomic radii after substitution. It is known that is 1.72 Å, while is only 0.79 Å.

Table 2.

The optimized results for occupation sites, cell parameters, volume, and total energy of Fe1.75H0.5SiO4.

.

To explore the relation among the mass density, lattice volume, and the hydrostatic pressure, the trends of the mass density, lattice volume as a function of the pressure are drawn in Fig. 2, respectively. As can be seen from the curve, the greater is the pressure, the smaller the volumes of the system Fe2SiO4 and Fe1.75H0.5SiO4 are, while the mass densities have the opposite trend. Compared with Fe2SiO4, the mass density of Fe1.75H0.5SiO4 slightly decreases, which apparently results from the reduction of total mass after substitution. In addition, the lattice volume of Fe1.75H0.5SiO4 reduces under the same pressure, which results from the decrease of atom radius after substitution. The greater is the pressure, the smaller the difference of the system densities Fe2SiO4 and Fe1.75H0.5SiO4 (6.7% at 0 GPa, 3.7% at 30 GPa) are.

Fig. 2. (color online) The trends of the density, cell volume as a function of the pressure for Fe2SiO4 and Fe1.75H0.5SiO4.
3.2. Elastic moduli of anhydrous and hydrous fayalite

Fayalite has an orthorhombic lattice with Pbnm space group, which has nine independent elastic constants. Bulk and shear moduli (K and G) are two important parameters of elasticity for fayalite, which can be obtained from elastic constants cij employing the following expressions:[58] In addition, the P and S wave velocities ( and ) of the system can be obtained by means of the elastic moduli K and G as follows: where ρ is the mass density of the system, which can be calculated by total mass and volume of the fayalite.

To discuss the effect of hydration on the elasticity of fayalite, the elastic constants of Fe2SiO4 and Fe1.75H0.5SiO4 at hydrostatic pressures from 0 GPa to 30 GPa are given in Table 3 together with the previous measured results. Firstly, our calculated elastic constants of Fe2SiO4 at 0 GPa are contrastively analyzed with the measured results under ambient conditions. The zero-pressure values are basically in agreement with the experimental data. The discrepancy with experimental results may be due to some temperature dependence. The present results were obtained under the condition of 0 K. In order to compensate for the influence of temperature, the elastic constants of Fe2SiO4 under 3 GPa pressure condition were obtained as shown in Table 3. It is obvious that the 3 GPa pressure correction brings a good agreement between our calculation and previous experiment at room temperature.[3133,55] As expected, the maximum absolute error is approximately 5 GPa for c11, the value deviates by less than 2%. The errors between the other elastic constants and the experimental data decrease to some extent.

Table 3.

Elastic constants of Fe2SiO4 and Fe1.75H0.5SiO4 at the different pressure conditions with the experimental data (in unit GPa).

.

In addition, all elastic constants of Fe2SiO4 gradually increase with increasing pressure. In order to make clear comparison, the longitudinal constant c33 of Fe2SiO4 and Fe1.75H0.5SiO4 as a function of pressure are plotted in Fig. 3 together with the measured results. As can be seen from Fig. 3, the elastic constant c33 is directly proportional to the pressure at least up to 30 GPa. Compared with the experimental data, the present data are slightly underestimated due to the nature of GGA. Unfortunately, the elastic constants of Fe1.75H0.5SiO4 have not been reported so far. Compared with Fe2SiO4, the elastic constants of Fe1.75H0.5SiO4 would slightly reduce under the different pressure conditions. The reason for the reduction may lie in the increase of the distance between O–O atoms which leads to the decrease of repulsive force between O–O atoms with the substitution of H atoms. In addition, as pressure increases from Fig. 3, the differences of c33 between Fe2SiO4 and Fe1.75H0.5SiO4 gradually lower. For Fe1.75H0.5SiO4, the difference is 54.08 GPa at 0 GPa, decreased by 23.85%, while the data (20.87 GPa) produces only 6.48% decrease at 30 GPa.

Fig. 3. (color online) The calculated and experimental elastic constant C33 of Fe2SiO4 and Fe1.75H0.5SiO4 under pressure.

The pressure dependence of elastic moduli of Fe2SiO4 and Fe1.75H0.5SiO4 are plotted in Fig. 4. For Fe2SiO4, the data of the bulk (K) and shear moduli (G) at 0 GPa were given in the previous experiment (K = 128 GPa–140 GPa and G = 49 GPa–55 GPa).[55] Our calculated moduli B and G are 129.15 GPa and 49.81 GPa, which is consistent with the measured data. It indicates from Fig. 4 that the bulk and shear moduli of Fe1.75H0.5SiO4 under different pressure conditions gradually reduce with the substitution of H atoms. At 0 GPa, the bulk/shear modulus of Fe1.75H0.5SiO4 is respectively 28.72% and 8.67% less than that of Fe2SiO4, while the bulk/shear modulus of Fe1.75H0.5SiO4 is 10.48% and 2.01% less than that of Fe2SiO4 when the pressure increases to 30 GPa. It indicates that as pressure increases, the differences of elastic modulus between Fe2SiO4 and Fe1.75H0.5SiO4 gradually lower. The reason is that the different degrees of volume change caused by pressure. It can also be seen from Fig. 4 that elastic modulus of Fe2SiO4 and Fe1.75H0.5SiO4 gradually increase with the increase of pressure up to 30 GPa.

Fig. 4. (color online) Pressure dependencies of bulk modulus (K) and shear modulus (G) of Fe2SiO4 and Fe1.75H0.5SiO4.
3.3. Effect of hydration on wave velocities of fayalite and its geophysical implication

Employing the mass density (ρ) and bulk and shear moduli (K and G), the wave velocities of hydrous fayalite under the different depths were calculated and analyzed employing Eqs. (3)–(4). As we all know, the velocity of seismic wave is very helpful to detect the composition and dynamics of the mantle. Figure 5 provides the compressional ( ) and shear ( ) velocities of Fe2SiO4 and Fe1.75H0.5SiO4 (fayalite with 2.3 wt% H2O). Compared with Fe2SiO4, the P ( ) and S ( ) velocities of Fe1.75H0.5SiO4 at 0 GPa is slower respectively by 8.5% and 1.1%, while the differences at 25 GPa reduce to 3.4% and 0.3%. It indicates that the differences of wave velocities between Fe2SiO4 and Fe1.75H0.5SiO4 are inversely proportional to the pressure. Owing to the larger pressure derivatives of K and G, the change trend of the pressure induced sound velocities for Fe1.75H0.5SiO4 is more obvious than that of Fe2SiO4.

Fig. 5. (color online) Calculated acoustic velocities ( and ) of Fe2SiO4 and Fe1.75H0.5SiO4 with depth.

Based on the seismic exploration data, a low velocity zone (LVZ) exists in the mantle region.[5962] However, the origin of LVZ has been debated so far. Several mechanisms have been applied to explain the origin of LVZ.[6365] The research suggests that the main reason for LVZ is the presence of water in mantle rocks. According to our data, adding 2.3 wt% water into fayalite results in the decrease of 3.4%–7.5% ( ) and 0.3%–3.4% ( ) from 0 GPa to 25 GPa, respectively. Our results are basically in agreement with the 2%–5% reductions of sound velocity obtained by the experimental observation in the LVZ. The minor difference between our data and experimental results, may be caused by ignoring the temperature factor during the our calculations.

3.4. Electronic structures of anhydrous and hydrous fayalite

To discuss the effect of hydration on the electronic properties of fayalite, firstly, the total and partial densities of states (TDOS and PDOS) are shown in Fig. 6 for Fe2SiO4 and Fe1.75H0.5SiO4 at 20 GPa. and the zero point of energy in every figure represents the Fermi energy . It is obvious that the TDOS at is 0 states/eV-f.u in Fig. 6. It indicates that hydrous fayalite still keeps the insulating characteristics, which agrees well with the experimental evidence of insulator.[66]

Fig. 6. (color online) The total densities of states (DOS) and the local contributions of particular atoms of Fe2SiO4 (a) and Fe1.75H0.5SiO4 (b).

For Fe2SiO4, the valence band is mainly dominated by the Fe-3d and O-2p states, whose width is approximately 9.1 eV. The width of conduction band is about 7.2 eV, which is dominated by the Fe-3d spin-down states. For Fe1.75H0.5SiO4, the valence and conduction bands are still provided generally by the Fe-3d and O-2p states, which are partially contributed by H-1s states. As is known to all, the more delocalization of electrons and the more width of valance bands, the stronger are interactions among atoms of the system. Compared with Fe2SiO4, the valence band of Fe1.75H0.5SiO4 is broader by about 0.7 eV, which indicate that the interactions among atoms could be strengthened with the substitution of Fe atom by H atoms. In addition, the large-scale overlap between Fe-3d spin-up states and O-2p states also can be observed from Fig. 6, which indicates that the electrons between Fe-3d and O-2p would be hybridized and bonded in the almost whole energy range of valence band. The top of valence band near originates from the Fe1-3d and the Fe2-3d spin-down orbits. Furthermore, the hybridization between O-2p and H-1s orbits mainly occurs at the bottom of valence band, which would probably enhance the stability of this system.

It is meaningful to study the dependence of the electric resistivity of hydrous fayalite and the external pressure. The band gaps of Fe2SiO4 and Fe1.75H0.5SiO4 are plotted in Fig. 7 with increasing pressure. It can be clearly observed from Fig. 7 that anhydrous and hydrous fayalite keep the insulation characteristics under the pressures up to 30 GPa, which is in agreement with the measured results.[56,67] For anhydrous fayalite, the band gap (2.21 eV) at 0 GPa is close to the previous calculated value 2.52 eV,[51] but our calculated gap is underestimated in comparison with the experimental value 4.22 eV obtained by spectroscopic technique.[68,69] In addition, the band gap of anhydrous fayalite gradually decreases with increasing pressure, which agrees well with the experiment conclusion.[19,56,66] For hydrous fayalite, the band gap width does not decrease, but somewhat increases under the same pressure. Our calculated data have the same trend with the theoretical results of hydrous forsterite,[44] but it seems contrary to the previous experimental conclusions.[19,70] It indicates that H atoms may enter into the interstitial sites of fayalite in the form of lattice defects rather than substituting Fe atom. In other words, the experimental conclusions may result from other physical mechanisms including occupying the interstitial sites for H atoms or cation enhanced conductivity.[44]

Fig. 7. (color online) The effect of pressure on band gaps of Fe2SiO4 and Fe1.75H0.5SiO4.

The calculated spin magnetic moments for anhydrous and hydrous fayalite are plotted in Fig. 8. For Fe2SiO4, the spin moment of Fe1 atom is , while that of of Fe2 atom is , which are less than the experimental values measured at 10 K ( for Fe1 and for Fe2).[71] Our data are closer to the previous theoretical predictions.[4951] The magnetic structure of Fe2SiO4 is complex, the present and previous simulations are only simplified as a collinear magnetic structure. The spin magnetic moments of the different Fe sites for anhydrous fayalite gradually decrease with increasing pressure, which has the same trend with the previous results.[50,51] Compared with Fe2SiO4, the magnetic moment of Fe1 atom for Fe1.75H0.5SiO4 slightly decreases, while the data of Fe2 atom basically keep constant. The spin moments of Fe1 atom and Fe2 atom of Fe1.75H0.5SiO4 linearly decrease with increasing pressure.

Fig. 8. (color online) The calculated spin magnetic moments of Fe2SiO4 and Fe1.75H0.5SiO4.
4. Conclusions

The stable structures, elasticities and electronic properties of Fe2SiO4 and Fe1.75H0.5SiO4 have been studied in detail by density functional theory. According to modeling and optimizing, two H atoms prefer to substitute Fe atom of Fe1 sites. As the pressure increases, the lowering of cell volumes and the increasing of mass densities for Fe2SiO4 and Fe1.75H0.5SiO4 are consistent with the experimental data. In addition, nine elastic constants and elastic moduli are proportional to the pressure for Fe2SiO4 and Fe1.75H0.5SiO4.

According to our calculations, adding 2.3 wt% water into fayalite (Fe1.75H0.5SiO4) reduces the wave velocities by 3.4%–7.5% ( ) and 0.3%–3.4% ( ) from 0 GPa to 25 GPa, which are basically in agreement with the 2%–5% reductions of sound velocity obtained by the experimental measurement in the LVZ. The present research suggests that the main reason for LVZ is the existence of water in mantle rocks. Based on the data of electronic structure, anhydrous and hydrous fayalite still keep the insulation characteristics under the pressures up to 30 GPa, which seems contrary to the previous experimental conclusions.

Reference
[1] Martin R F Donnay G 1972 Am. Mineral. 57 554
[2] Beran A Putnis A 1983 Phys. Chem. Miner. 9 57
[3] Tarits P Hautot S Perrier F 2004 Geophys. Res. Lett. 31 265
[4] Rosa A D Mezouar M Garbarino G Bouvier P Ghosh S Rohrbach A Sanchez-Valle C 2013 J. Geophys. Res. 118 6124
[5] Balan E Blanchard M Lazzeri M Ingrin 2014 Phys. Chem. Miner. 19 460
[6] Bindi L Nidhi M Tsuchiya J Irifune T 2014 Am. Mineral. 99 1802
[7] Ghosh S Schmidt M W 2014 Geochim. Cosmochim. Acta 145 72
[8] Ragozin A L Karimova A A Litasov K D Zedgenizov D A Shatsky V S 2014 Russ. Geol. Geophys. 55 428
[9] Ohira I Ohtani E Sakai T Miyahara M Hirao N Ohishi Y Nishijima M 2014 Earth Planet. Sci. Lett. 401 12
[10] Xu Z Zheng Y F Zhao Z F Gong B 2014 Geochim. Cosmochim. Acta 143 285
[11] Ohtani E 2015 Chem. Geol. 418 6
[12] Green H W II Burnley P C 1989 Nature 341 733
[13] Frost D J Dolejš D 2007 Earth Planet. Sci. Lett. 256 182
[14] Cornwell D G Hetényi G Blanchard T D 2011 Geophys. Res. Lett. 38 239
[15] Schmandt B Jacobsen S D Thorsten W Becker T W Liu Z X Dueker K G 2014 Science 344 1265
[16] Karato S 2013 Phys. Earth Planet. Inter. 219 49
[17] Dai L D Karato S I 2014 Phys. Earth Planet. Inter. 237 73
[18] Karato S 1990 Nature 347 272
[19] Wang D J Mookherjee M Xu Y S Karato S I 2006 Nature 443 977
[20] Huang X Xu Y S Karato S 2005 Nature 434 746
[21] Manthilake M Matsuzaki T Yoshino T Yamashita S Ito E Katsura T 2009 Phys. Earth Planet. Inter. 174 10
[22] Yoshino T 2010 Surv. Geophys. 31 163
[23] Pearson D G Brenker F E Nestola F McNeill J Nasdala L Hutchison M T Matveev S Mather K Silversmit G Schmitz S Vekemans B Vincze L 2014 Nature 507 221
[24] Jacobsen S D Jiang F Mao Z Duffy T S Smyth J R Holl C M Frost D 2008 J. Geophys. Res. Lett. 35 L14303
[25] Wang D J Sinogeikin S V Inoue T Bass J D 2006 Geophys. Res. Lett. 33 L14308
[26] Mao Z Jacobsen S D Jiang F M Smyth J R Holl C M Frost D J Duffy T S 2008 Earth Planet. Sci. Lett. 268 540
[27] Kudoh Y Takeda H 1986 Physica B + C 139–140 333
[28] Andrault D Bouhifd M A Itié J P Richet P 1995 Phys. Chem. Miner. 22 99
[29] Zhang L 1998 Phys. Chem. Miner. 25 308
[30] Fukizawa A Kinoshita H 1982 J. Phys. Earth 30 245
[31] Sumino Y 1979 J. Phys. Earth 27 209
[32] Isaak D G Graham E K Bass J D Wang H 1993 Pure Appl. Geophys. 141 393
[33] Graham E K Schwab J A Sopkin S M Takei H 1988 Phys. Chem. Miner. 16 186
[34] Isaak D G Anderson O L Goto T Suzuki I 1989 J. Geophys. Res. 94 5895
[35] Liu Q Liu W Whitaker M L Wang L P Li B S 2010 Am. Mineral. 95 1000
[36] Piekarz P Jochym P T Parlinski K Lażewski J 2002 J. Chem. Phys. 117 3340
[37] Gillan M J Alfe D Brodholt J Vocadlo L Price G D 2006 Rep. Prog. Phys. 69 2365
[38] Liu L Du J G Zhao J J Liu H Wu D Zhao F L 2008 Comput. Phys. Commun. 179 417
[39] Yu Y G Vinograd V L Winkler B Wentzcovitch R M 2013 Phys. Earth Planet. Inter. 217 36
[40] Wang K Brodholt B Lu X C 2015 Geochim. Cosmochim. Acta 156 145
[41] Zhang X L Wu Y Y Shao X H Lu Y Zhang P 2016 Chin. Phys. 25 057102
[42] Gu J B Wang C J Zhang W X Sun B Liu G Q Liu D D Yang X D 2016 Chin. Phys. 25 126103
[43] Zhang H J Li S N Zheng J J Li W D Wang B T 2017 Chin. Phys. 26 066104
[44] Liu L Du J G Zhao J J Liu H Gao H L Chen Y X 2009 Phys. Earth Planet Inter. 176 89
[45] Tsuchiya J 2013 Geophys. Res. Lett. 40 4570
[46] Tsuchiya J Tsuchiya T 2009 J. Geophys. Res. 114 B02206
[47] Panero W R 2010 J. Geophys. Res. 115 B03203
[48] Wang D J Karato S I Liu Z Y 2012 Geophys. Res. Lett. 39 L06306
[49] Cococcioni M Dal Corso A de Gironcoli S 2003 Phys. Rev. 67 094106
[50] Jiang X Guo G Y 2004 Phys. Rev. 69 155108
[51] Stackhouse S Stixrude L Karki B B 2010 Earth Planet. Sci. Lett. 289 449
[52] Liu L Du J G Liu W Zhao J J Liu H 2010 J. Phys. Chem. Solids 71 1094
[53] Kresse G Furthmüller A 1996 Phys. Rev. B 54 11169
[54] Cococcioni M de Gironcoli S 2005 Phys. Rev. 71 035105
[55] Speziale S Duffy T S Angel R J 2004 J. Geophys. Res. 109 B12202
[56] Williams Q Knittle E Reichlin R Martin S Jeanloz R 1990 J. Geophys. Res. 95 21549
[57] Kudou Y Kuribayashi T Kagi H Inoue T 2006 J. Miner. Petrol. Sci. 101 265
[58] Hill R 1952 Proc. Phys. Soc. 65 349
[59] Hill D P 1972 Geol. Soc. Am. Bull. 83 1639
[60] Grand S P Helmberger D V 1984 Geophys. J. R. Astronom. Soc. 76 399
[61] Priestley K F Cipar J Egorkin A Pavlenkova N I 1994 Geophys. J. Int. 118 369
[62] Thybo H Zhou S Perchuc E 2000 Geophys. Res. Lett. 27 3953
[63] Artemieva I M 2003 Earth Planet. Sci. Lett. 213 431
[64] Gung Y Panning M Romanowicz B 2003 Nature 422 707
[65] Mierdel K Keppler H Smyth J R Langenhorst F 2007 Science 315 364
[66] Shankland T J 1975 Phys. Earth Planet Inter. 10 209
[67] Mao H K Bell P M 1972 Science 176 403
[68] Smith H G Langer K 1982 Am. Mineral. 67 343
[69] Burns R G 1970 Am. Mineral. 55 1608
[70] Xu Y S Poe B T Shankland T J Rubie D C 1998 Science 280 1415
[71] Fuess H Ballet O Lottermoser W Ghose S Coey J M D Salje E 1988 Arthrosonographie Berlin Springer